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Abstract 



A development of an inverse first-order divided difference operator for functions of several variables is presented. 

■ Two generalized derivative-free algorithms builded up from Ostrowski's method for solving systems of nonlinear 
equations are written and analyzed. A direct computation of the local order of convergence for these variants of 
Ostrowski's method is given. In order to preserve the local order of convergence, any divided difference operator 

■ is not valid. Two counterexamples of computation of a classical divided difference operator without preserving 
' the order are presented. A new divided difference operator solving this problem is proposed. Furthermore, a 

computation that approximates the order of convergence is generated for the examples and it confirms in a 
numerical way that the order of the methods is well deduced. 
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^ ■ 1 Introduction 

' In this paper, we consider iterative methods to find a simple root of a system of nonhnear equations 

F{x) = 0, (1) 

where F : D C R™ — > R"" and D is an open convex domain in R™. We assume that the solution of ([T]) is 
a G I? at which F'{a) is nonsingular. 



The goal of this paper is to give a generalization for several variables of a given Newton's type iterative method using 
' an approximation of the derivative by means of a central difference and preserving the local order of convergence. We 
center our analysis on Ostrowski's method ^ and the improvement but similar analysis should be performed 
on other methods. 

We would like to mention that the matrix of the divided difference operator, used generally in practice [J [21 [21 IH 
[51 [5] , does not preserve the local order of convergence in general. We present two counterexamples of this fact and 
we propose a new divided difference operator preserving the local order. 

After setting our iterative method we recall Ostrowski's iterative function [7] for a single nonlinear equation, 
f{x) — 0, with fourth-order that is given by 

/(^) 

y ^ ^- TTTT' 
/ W 

X = y-fif{y), (2) 

where / : R — > R, A* = 1/(2 [y,x;f] — f'{x))) and [y,x;f] — {f{y) — f{x))/{y — x) is the divided difference of 
/ at the points x and y. 

A variant of Ostrowski's method with sixth-order convergence is proposed in [S^. This improvement has the 
following scheme 
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/(^) 

/ W 

z = y-t^f{y), (3) 

X = z-fif{z). 

This method improves the local order of convergence of Ostrowski's method with an additional evaluation of the 
function. Other works related with the improvement of the local order of convergence of ([2]) can be found in 

In [13] a well-known technique (see [TS]) that consists in replacing the first derivative of ([5]) and ([3]) by a central differ- 
ence is used in order to obtain a new method that preserves the fourth and the sixth convergence order respectively. 
Notice that these iterative methods arc derivative free, since f'{x) is approximated by [a; -I- f{x),x — f{x);f]. 

The paper is organized as follows: We start in Section 2 with a generalization for several variables of the Ostrowski's 
method ^ and the improvement ^ using an approximation of the derivative by means of a central difference and 
preserving the local order of convergence. Section 3 is devoted to the computational efficiency index performing 
a comparison for the three iterative methods. By means of two counterexamples, in Section 4, we emphasize the 
necessity of a new way to compute the matrix of divided difference operator since in many examples the theoretical 
one and the classical way used until now to calculate it do not coincide. We finalize the paper with some concluding 
remarks in Section 5. 



2 Local order of convergence 

We assume that F has, at least, fourth-order derivatives with continuity on D. We consider the first divided 
difference operator of F on R™ as a mapping [■, ■ ; F] : D x D C R™ x R"" — > (R™) which is defined by (see 
[l6l m [17] and the references therein) 

[x + h,x;F]= I F'{x + th)dt, V(a;, /i) G R"" x R"". (4) 

Developing F'(x + th) in Taylor's series at the point x and integrating we obtain 

[x + h,x;F] = F'{x) + i F"{x) F"'{x) + 0{h^). (5) 

Taking into account that e = x — a we develop F{x) and its derivatives in a neighborhood of a. Namely, assuming 
that [F'{a)]~ exists, we have 

F{x)^r{e + A2e^+ Ase^+Oie^)), (6) 

where T = F'(a) and Ap ^ ^ T-^F^P^a) G ^p(R™,R'"), p = 2, 3. From (O the derivatives of F(x) can be 
written as 

F'{x) = r (/ + 2 ^2 e + 3 A3 + 0(6^)) , F"{x) = T (2 A2 + 6 A3 e + 0{e^)) , F"'{x) = T {6 A3 + 0(e)) . (7) 

Setting y ~ X + h and E ^ y ~ a, we have h = E — e. Replacing the previous expressions (O into ([5]) we get 

[y,x-F] = F'{a){l+ A2{E + e)+ A3{E^+Ee + e')+ O^iE.e)) , (8) 

where 03{E, e) represents the terms of order higher than three. We say that a function depending on E and e is 
an OsiE, e) if it is an OiEi" e«i) with go + gi = 3, qt>0, i = 0, 1. 

In all our analysis we must consider the operator [x + F(x), x — F(x); F]. From x + F{x) — a ^ e + F{x) and 
X — F{x) — a = e — F{x), if we replace in 1^ E hy e + F{x), e by e — F{x), and F{x) by then we will get 

[ X + F{x), X - F{x); F] = F (/ + 2 A2 e + A3 (3 + e 2) + 0(6^)) , (9) 

where e = F'{a) e. Expanding in formal power developments the inverse of the preceding operator ([9]) in terms of 
e and e we get 

[x + F{x),x ^ F{x);Fy^ = {I - 2 A2 e ^ Asi^e^ + e^)+ AAje^ + 0{e^)) T'K (10) 



The first step of the variant of Ostrowski's method is Newton's method using a central divided difference. That is, 
we have the following iteration function 

y = $o(a;) = x~ [x + F{x),x ~ F{x)]FY^ F{x). (11) 

Subtracting the root a from both sides of pT|) and considering ([5]) and (|10p . the corresponding vectorial error 
equation is given by 

E = y-a = e- {I -2A2e + 0{e^)) {e + A2e^ + 0{e^)) 

The second step of the variant of Ostrowski's method is defined by the following iterative function 

z ^ <^>i{x,y) = y- i^F{y), (12) 

where ly = {2 [y,x;F] — [x + F{x), x ~ F{x); F])~^ . Note that from the developments given in (|8]) and ^ we 
have 

jy-^ = T{I + 2A2E+ A3{2E^ + 2eE~e'^ -e^)+Oie^)) 
- r{l + 2Ale^- A3{e^ + e^)+0{e^)), 

and applying again the computation of the inverse operator we have in this case 

V = {I -2Ale^ + A^ie^ + e^)+0{e'')) T-\ 

Subtracting the root a from both sides of (|T2|) the vectorial error equation for the iterative method defined in ((TT|) 
and ([HI) is 

e^z-a = E- {l-2Ale^+ A3{e^ + e^)+0{e^)){E + A2E^+ 0{E^)) 

= {Al- A3A2) e^- A3e^A2e^+ Oie'). (13) 

Summarizing the preceding results we are in a condition to state the following theorem. 

Theorem 2.1 The iterative method given by the iteration function $1 defined in Ul\) and US\) is at least of fourth 
local convergence order and its vectorial error difference equation is given in il3\) . 

Moreover, if we introduce a new step with an additions evaluation of the function F in the following way 

X = ^2{x,y,z) = z- iyF{z), (14) 

then the error becomes 

X-a ^ e- {I -2Ale^ + A3{e^ + e^) + 0{e^)) {e+ Oie"^)) 

= [(2 Aj - A3) - A2 e 2] [{Al - A3 A2) e' - A3 e' A2 e'] + 0{e'), (15) 

and the local error of convergence of the iterative method defined by the steps (|lip . ([T^ and ([H)) is at least of 
sixth order. We present the preceding results in the following theorem. 

Theorem 2.2 The iterative method given by the iteration function $2 defined by ill]) . il2\) and ^14\ l is at least of 
sixth local order of convergence and the vectorial error difference equation is il5\) . 



3 Computational efficiency and comparison between the methods $i 

The computational efficiency index {CEI) and the computational cost per iteration (C) are defined as 

CEI{^i,m,l) = , and C{fi,m, £) a{m) + p{m, £), (16) 

where a{m) represents the number of evaluations of the scalar functions (Fi, . . . , Fm) used in the evaluation of F, 
[x — F{x),x + F{x);F] and [y,x;F]. The number of products needed per iteration is represented by p{m,£). In 
order to express in terms of products the value of C(/i, m, £) a ratio /i > between products and scalar functions 
evaluations and a ratio £ > 1 between products and quotients are required (see [T71 [H]). A first and classical 
definition of difference divided operator (|121 HI HOI 131 HI]) used in the computations is 

Definition 3.1 

[y^x'^F]^^] = {Fiiyi, . ■ . ,yj-i,yj,Xj+i, . . . ,Xr,r) - F,{yi, . . . ,yj_i,Xj,Xj+i, . . . ,Xr,r)) /{yj-x^), 1 < i,j < m. (17) 



This operator defined is a divided difference of F on the points x and y because it is a bounded linear operator 
which satisfies the condition [331 US] 

[y,x-F](y~x)^F{y)-F{x). (18) 



When we evaluate F in any iterative function we calculate m component functions and if we compute a divided 
difference then we evaluate m{m — l) scalar functions, where F{x) and F{y) are computed separately. We must add 
quotients from any divided difference. In order to compute an inverse linear operator we have m(m— 1) (2m— 1) / 6 
products and TO(m — l)/2 quotients in the decomposition LU and m{m — 1) products and m quotients in the 
resolution of two triangular linear systems. Taking the previous considerations, for the first iterative function $o 
we have a{m) = m? + 2m, p{m) = m(2m^ + 3m — 5)/6 + £m{3m + l)/2 we obtain 

Co = m(m + 2)/i+ — (2m2 + 3m- 5) + — (3m + 1) and CEIq ^ 2^/^^° . (19) 
6 2 

In an analogous way for $i we get 

Ci = 2m(m+l)/i+^(2m2 + 3m-5) + £m(3m + l) and CEh = i^^'^^ (20) 

o 

Finally, for the third iterative method we take a{m) — m{2m + 3), p{m) — m{2m? + 6to — 8)/3 + £m{3m + 2) and 
the cost and the computational efficiency index are 

C2 = m(2m + 3)Ai+— (2m2 + 6TO-8)+^m(3m + 2) and CEl2 = 6^^'^^. (21) 
3 

In order to compare the efficiency of the iterative methods $i, < z < 2, we define the ratio 

^ ^ log CEF{fi,m,i) ^ log(pj) Cj{^i,m,i) 
''^ log CEIj {^i,m,e) \og{pj )Ci{^i,m,i)' 

For Rij > 1 the iterative method $i is more efficient than Notice that when if Ri j = 1, then we have the 
border between the two computational efficiencies. We study this equation for each pair of methods to obtain the 
following theorem: 

Theorem 3.1 For all £ > 1 and m>2 we have: 

CEh > CEh > CEIq 

Proof. 

1. In the comparison {i,j) = (2, 1) we have 

^ log p2 2m^ + 6m/i + 3m + 9i?m + 6/i + 3f — 5 
^ ^ log pi 2m? + 6m/i + 6m + 9^m + 9/i + 6i? — 8 

The equation R21 — 1 has one vertical asymptote for m = 0.7095 in the (m,/x) variables and for all t > 1. 
Furthermore, p is always negative for m > 2 and for all £ > 1. We conclude that CEI2 > CEIi since i?2 1 > 1 for 
m > 2, /i > and ^ > 1. 

2. In the comparison {i,j) — (1,0) we have 

^ _ 2m2 + Qmp + 3m + Mm + 12^ + 3^-5 
^ ^ 2m2 + Qmp + 3m + 9^m + 6/^ + 3^ — 5 

Equation i?io = 1 is equivalent to 6/j, = for all ^ > 1 and m > 2. Finally, we conclude that CEIi > CEIq since 
i?io>lform>2, ^>0 and £>1. The proof is completed. □ 



3.1 Numerical analysis 

The numerical computations were performed on MPFR library of C++ multi-precision arithmetics [53] with 4096 
digits of mantissa. All programs was compiled by gcc(4.3.3) for i486-linux-gnu with libgmp (v. 4. 2. 4) and 
libmpfr (v. 2. 4.0) libraries on an Intel(8)Xeon E5420, 2.5GIIz and 6MB cache. For this hardware and software 
the computational cost of the quotient respect to the product is £ — 2.5 (see Tabled]). Within each example the 
starting point is the same for all methods tested. The classical stopping criteria ||e/|| = ||a;/ — ajj > 0.5 • 10"*^ and 
||e/+i|| < 0.5 • 10"% with e = 4096, is replaced by 

Ei = ttJ^^ > 0.5 ■ 10"" and Ej+i < 0.5 • 10"", 



A 



where ei = xj — xj^i and 77 = e. Notice that this criterium is independent of the knowledge of the root (see 
[25]). Furthermore, in aU computations we have substituted the computational order of convergence (COC) [26] 
by an approximation (ACOC) denoted by p [53] and defined as follows 

ln£;/ 
^' \nEi_i ■ 

According to the definition of the computational cost (1161) . an estimation of the factor /i is claimed. To do this, we 
express the cost of the evaluation of the elementary functions in terms of products, which depends on the machine, 
the software and the used arithmetics fM] [27] . In Table [T] an estimation of the cost of the elementary functions in 
product units is shown, where running time of one product is measured in milliseconds. 

Table 1: Computational cost of elementary functions computed with a program write in C++, compiled by gcc(4.3.3) for 
i486-linux-gnu using libgmp (v. 4. 2. 4) and libmpfr (v. 2. 4.0) libraries on an Intel®Xeon E5420, 2.5GHz, 6MB cache, 
where x — \/3 — 1 and y = \/5. 

Arithmetics x*y x/y y/x exp{x) ln(a;) sin(a::) cos(x) arctan(a;) 

C++ MPFR 4096 digits 0.1039 ms 2.5 1.7 87.8 66 116 113 228 



We analyze three examples. The first one does not present any problem when we apply the definition 13.11 given 
in (|17|) . The second and third examples, presented in the next section, are in fact counterexamples. That is, they 
need the new definition 14.11 given in ([25]) to attain the local convergence order that theoretical results state. 

Tables [2]^ show the results obtained for the iterative methods < fc < 2. In each table we show the necessary 
iteration number /, the computational cost Ck in terms of products, the computational efficiency index CEIk 
([r9| -([2T ]) . the measure time factor TF defined as TF — l/lognjCEIk, the necessary time T in milliseconds to 
reach the /th-iteration, and the computed value of ACOC with its higher bound Apfe, where the local order is 
computed by pfe = ± Apk |25j . Furthermore, the last column shows the number of correct decimals in Xj, say q. 

We finalize this section with a first numerical comparative study in the following system of nonlinear equations we 
consider the system F{x) = with five exponential equations defined by 





X2 ' 


\^X3 + 


2:4 H 


- - 








Xi - 


l-a;3 + 




- - 






F{xi,X2,X3,X4,X5) = < 


Xl - 


^X2 + 


X4 H 


- 2^5 - 




(22) 




Xi - 


t-xz + 


2:3 H 


- Xs - 










t-xa + 




- a;4 — 






The results shown in table [2] where xq 




-2.1,- 


-2.1,6.4,6.4, 


—2.1)* tends towards the 



where {m, fi) = (5, 87.8). 

root a « (-2.153967996 . . . , -2.153967996 . . . , 6.463463374 . . . , 6.463463374 . . . , -2.153967996 ...)*. We verify the 
assertion of the theorem 13.11 that is, the efficiency of iterative method $2 is better than the efficiency of $1 and 
$1 one is better than $0 one. Moreover, we obtain the local order of convergence given in theorems 12. II and 12.21 



Table 2: Numerical results for the non linear system 
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CEI 


TF 


T 


/5± Ap 




1 




11 


3223.0 


1.000215086 


10706.57 


3730 


2 + 3.04-10- 


24 


3493 


*1 


5 


5568.0 


1.000249006 


9248.26 


2940 


4 + 7.89-10- 


11 


1112 


$2 


4 


6039.5 


1.000296717 


7761.36 


2540 


6 + 2.00-10- 


-7 


1191 



4 First-order divided difference operators preserving the local order 

This section is devoted to the preservation of the local order of convergence, of a given method, when we replace 
the derivatives by divided differences. The most usual definition in practice of divided difference is given in (jl7p . 
However, we will need a new definition of divided differences instead of definition 13.11 given in ([T7| to get the local 
order wanted when we apply algorithms fc = 0, 1, 2. We present two counterexamples of this fact. 
Case 1. The first numerical counterexample presented consists in the following 2x2 system of nonlinear equations 



f x2 + a;2 - 9 = 0, 
F(xi,X2) = \ (23) 

[ Xi X2 - 1 = 0, 

where we have ^ = 1.5. When we test the behavior of the family of iterative methods < /c < 2, with initial 
value xq = (3.0,0.4)* towards the root a w (2.98118805,0.335436739)*, the computation of ACOC give 2, 3 and 4 
respectively, against the result wanted. Furthermore, a comparison of the expression (j4]), taking into account the 
definition 13. 11 gives the following result 



/•I , / 2x1 + hi 2x2 + h2 \ 

/ F'{x + th)dt^\ \=^\x + h,x;F 



(-^j / 2x1 + hi 2x2 + h2 
X2 + h2/2 xi+hi/2 J ' ' ■ ■ ' 1 2;i + /ii 



where x = (xi, X2)'^, /i = (/ii, /12)"'" and t G R. Due to Potra [BEH] we have the following necessary and sufficient 
condition to characterize the divided difference operator by means of a Riemann integral: 

Theorem 4.1 If F satisfies the following Lipschitz condition \\[x,y;F] — [u,v]F]\\ < H (\\x — u\\ + \\y — v\\), then 
equality ^ holds for every pair of distinct points (x + /i,x) G D x D if and only if for all {u,v) E D x D with 
u ^ V and 2v — u E D it is satisfied the relation 

[u,v;F] = 2 [u, 2v - u; F] - [v,2v - u; F] . (24) 



It can be checked that the function considered in (1^^ does not hold (IMl) . We need a new definition of divided 
differences instead of definition 13.11 given in ([T7| to get the local order wanted when we apply algorithms k = 
0, 1, 2, in this case. We take the following way to compute the divided difference operator: 

Definition 4.1 

[y,x;F]f^^{F,{yi, . . . , yj_i, y^, Xj+i, . . . ,x™) - Fi{yi,.. . , ?/j_i, Xj, Xj+i, . . . ,x„) 

+Fi(xi, . . . ,Xj-i,yj,yj+i, . . . ,ym) - F,{xi, . . . ,Xj-i,Xj,yj+i, . . .,ym))/ (2 {yj ~ x/)) , i<»j<m. (25) 



Notice that this operator defined in ([25]) is not only a divided difference of F on the points x and y because it 
accomplishes the condition (jl8p but also is a symmetric operator: [y,x;F] = [x^y]F]. If we use the definition 14.1 
we have to double the number of evaluations of the scalar functions in the computation of [y, x; F]. 

A development of Taylor of the expressions considered in the previous definitions gives 

l^l ^ ni ^ m 

/ DjF,{x + th)dt = DjF,{x) + -y^DkjF,{x)hk + -y^DkijF,{x)hkhe + 0{h^), 
•^0 ^ k=i ^ 

[x + h,x; F]^ = DjF,{x) + ^ Dfe,F,(x) hk + -D,,F,{x) h, + Oih''), 

k=l 

[x + h,x;F]lf = DjF,{x) + -Y,DkjF,{x)hk + - ^ DkijF,{x)hkhi + -D,j,F,{x)h] 

k=l k,t=l 

k + l<2j 

_^ m 

+ 1 Dk:ejF,{x)hkhe + 0{h^^ 

k + l>2j 

From the previous expansions we obtain 

I DjF,{x + th)dt-[x + h,x]F]'lf ^ 0{h), and / DjF,{x + th) dt - [x + h,x; F]'lf = Oih"^), 
Jo Jo 

and the definition 14.11 has higher precision. Moreover, if we use the computational definition 13.11 given in ([T7)) . then 
the operator = 2 [y, x; F] — [x + F{x), x — F{x); F], can be written as 



=T[2I+ 2B2e + Oie/)j - T[^I + 2^2 e + 0(e") j , 

where i?2, -B2 G ^2 (R-™, R-™): the set of bounded 2-linear operators, but B2 ^ A2 and B2 ^ A2, since in (|23p we 
have DkjFi{x) ^ 0. This fact produces a loose of convergence order in many algorithms. Notice that in this case 



it appears the term in e and the increase of local order of iterative method $i respect to $o (HB) is only one. 
On the other hand, if we consider the computational definition 14.11 given in (j25p , when we develop the operator 
then the difference with the theoretical result appears in the terms of second degree in e, while in first order 
are identical. Namely, 



ly-^ =r(2I+ 2^26+ 2{B3 + Al) + 0{e^)'^ - r(/+ 2^26+ B3 {3e^ + + Oie^)^ 



(26) 



where S3, B3 e Jfa (R™, R™). 

Using definition 14.11 for computing the divided difference operator we have to add more evaluations of functions. 

('2) ('2) 

The new CEI (C) for $fc indicates by CEIl ' {Cl '), fc = 1, 2. Note that for $0 it is not necessary the new definition. 



(2) 



4mV+ — (2m2 + 3m- 5) +£m(3m+ 1), CEl[^'' = 3'^/'^^ , CEl[^^ = i^^'^^ 



Cf^ = m(4TO + l)^ + — (2m2 + 6m-8)+^m(3m + 2), Ci^^^^ = 4 ^/^^^'^ CEI^/^=6^/^- 



(1) 



r(2) 



-.(2) 



With these new computational efficiencies we obtain the following theorem: 
Theorem 4.2 For all m > 2 and £ > 1 we have 



but CEIq can be greater than CEl'^^ and CEl2^' . That is 

• Form = 2 and for all £ > 1 we have CEIq = CEl[^^ and CEI^ > CEI^^^ . 

• For all m>2 and£>l we have CEIq > CEl[^\ 

• For other comparisons it depends on the values of 11, m and £ (see proof). 



CElf^ > CEI^^\ 



50 - 
45 - 
40 - 
35 - 



Figure 1: Boundary curve = G^q^' for £ = 2.5 (f27)) 



Proof. 

1. In the comparison {i,j) = (2, 1) we have 



^(2 2) _ log P2 2m^ + I2mn + 3m + Urn + 3f - 5 

~ log pi 2m? + 12to^ + 6m + Um + 3^ + 6£ - 8 ' 



f2 21 

The equation i?2i ^ 1 ^as one vertical asymptote for m = 0.8548 in the (m,/i) variables and for all £ > I. 
Furthermore, is always negative for m > 2 and for all £ > 1. We conclude that CEI2 > CEI{ since i?2i > 1 
for m > 2, ^ > and £>l. 

2. In the comparison {i,j) — (1,0) we have 



(2 1) _ 2m2 + 6m^ + 3m + 9£m + 12/i + 3^-5 



-^10 ^ 



2m2 + 12mAt + 3m + 9£m + 3£ - 5 



7 



(2 1) 

Equation R\q = 1 is equivalent to equation m = 2 for all £ > 1 and /i > 0. Finally, we conclude that CEIq > 
CEl[^^ since r'-^q^ < 1 for m > 2, ^ > and £ > 1 and CEIo = CEl[^^ for m = 2. 
3. In the comparison = (2, 0) we have 



^(2 1) _ log P2 2to2 + 6mfi + 3m + dim + 12/i + 3£ - 5 
^° ~ log po 4m2 + 24m/i + 12m + 18im + 6/x + 12£ - 16 ' 

(21) 

The particular boundary i? 2 1 expressed by /i written as a function of £ and m is 



,(2 1) _ 1 2qm^ + 3{3q£ - r)m - 3r£ - {2q - 3r) 
3 2rm — {7q + 3r) 



^(2 1) _ ^ ~<~ ■J\^q< — I /III — oif. — y^q — ^1 ) 

'-'2,0 — o KZZ , ' 



(2 1) 

where q = ln(3/2) and r — ln(8/3). The function G2 2 presents a vertical asymptote for m = {7q + 3r)/(2r) = 

f2) 

2.9468 . . . and it is always negative for m = 2 and £ > 1. In this case, we conclude that CEI^ > CEIq since 
R^^Q^ > 1 for TO = 2, /i > and ^ > 1. 

(2 1} 

Finally, for m > 2 the boundary /j, = q separates two regions which extend for increasing values of m and /i 
(see figure [1]). 

4. In the comparison = (2, 2) with different local order and cost we have 

^(2 1) _ log 6 2to2 + 6m^ + 6m + 9£m + 9/i + 6£ - 8 
log4 2m2 + 12m^ + 6m + 9^m + 3^ + 6£-8' 

(21) 

Here the boundary i?2 2 1 expressed by /i as a function of £ and m can be written as 



^2,2 -o — , ■ i^^J 



^(2 1) _ 1 2qm'^ + 3q{3£ + 2)m + 6q£ 
3 2rm - (5g + 2r) 

f2 1) 

This function, /i = 2 7 has a vertical asymptote for m = (5q + 2r)/(2r) — 2.0334 . . . and it is always negative 
for TO = 2 and £ > 1. In this case, we conclude that CEI^^ > CEI^^ since R^22' > 1 for m = 2, /i > and £>\. 
Finally, for m > 2 the boundary fj, = G2 2 separates two regions which extend for increasing values of m and /i 
(see figure [2]). 

5. In the comparison = (1, 1) with different local order and cost we have 

^(2 1) _ log 4 2to2 + 6m^ + 3m + 9^to + 6/x + 3^ - 5 
~ log3 2to2 + 12TO/Z + 3to + 9^TO + 3£ - 5 

(2 1) 

The particular boundary R\ = 1 expressed by p, written as a function of £ and to is 

(2 1) 1 2sTO2 + 3s(3£+l)m + 3s£-5s 

^1,1 — To 77 ^ 7 ' ^^^^ 

12 [t — s)m — t 

(2 1) 

where s — ln(4/3) and t = In 2. The function G\ ^ has a vertical asymptote for m = t/{t — s) = 1.7095 . . .. The 
boundary /i = G^^^^ separates two regions that are shown in figure [3l 
The proof is complete. □ 

Table 3: Numerical results for the non linear system (|23[) . 







C 


CEI 


TF 


T 






q 


$0 


11 


32.5 


1.021556664 


107.96 


40 


2 ±2.61 • 


10" 


-04 


3334 




7 


59.0 


1.018794991 


123.66 


40 


3 ±4.07 • 


10- 


-40 


2908 




5 


65.0 


1.021556664 


107.96 


30 


4 ±5.86- 


10 


-13 


1951 




5 


69.0 


1.020294410 


114.61 


35 


4± 1.70• 


10- 


-23 


1384 




4 


75.0 


1.024177781 


96.38 


20 


6± 1.90- 


10- 


-08 


2392 



Table [3] shows the same columns that table [2] but taking into account the two definitions of divided difference 
operator. Note that for this case 2, we have (TO,/i,£) = (2,1.5,2.5) and the theorem 14.21 is accomplished. That is, 

CEI^ > CEl[^^ = CEIo > CEl[^^ and CEI^^^ > Ce4^\ 



Q 



60 p 
45 - 



15- 
10 - 

5 - 

oL 



CEfl>>CEP> 



Figure 2: Boundary curve = G^^' for £ = 2.5 ([28] 
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Figure 3: Boundary curve fi = i^' for £ = 2.5 ([29] 



Observe that the values of ACOC in table [3] confirm the theoretical local orders of convergence with a high precision. 

Case 2. In this case we consider an example where the computational definition of the divided difference operator 
presents non zero values for and (see ([26])), opposite to the case 1. 

Xi — cos 2 a;^ — ^ = 0, I < i < m, (30) 

V j=i I 

where m = 3 and /.t = 113.3. We test the behavior of the iterative methods < fc < 2, say pT|) . (|12l) and 
(ITl) . with the initial value xq = (0.4,0.4,0.9)* towards the root a « (0.5438500415,0.5438500415,0.9957781534)*, 
considering the two definitions of the divided difference operator presented in this paper. 

Table HI shows the same columns that table [31 Notice that the values of ACOC confirm the theoretical local 
orders of convergence of the iterative methods presented with a high precision. In this case 3, we have (m,/i,€) = 
(3, 113.3,2.5) and the theorem 14.21 is accomplished. That is, 

CS/f ^ > CS/o > CS/f ^ > CS/f ^ and CEI^^'^ > CEI^^\ 



Table 4: Numerical results for the non linear system ((30]). 
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C 


CEI 


TF 


T 






q 




13 


1748.0 


1.000396616 


5806.73 


2330 


2 ± 1.58 • 


10" 


15 


2575 




8 


2816.2 


1.000390181 


5902.48 


2310 


3 ±9.56 • 


10- 


35 


2549 




6 


4175.8 


1.000332038 


6935.85 


2580 


4 ±5.88 • 


10- 


07 


2517 


$w 


6 


3169.6 


1.000437468 


5264.59 


2000 


4 ±5.35 • 


10- 


08 


1514 




4 


4529.2 


1.000395680 


5820.46 


1870 


6 ±7.52- 


10- 


06 


725 



5 Concluding remarks 

Ftom Ostrowski's method witch is a root finder for nonlinear equations we present two modified algorithms 
derivative-fi'ee for solving systems of nonlinear equations. A development of an inverse first-order divided dif- 
ference operator for functions of several variables is analyzed. Moreover, we study the computational efficiency 
index and a comparison for the three iterative methods is given. A similar analysis should be performed for other 
Newton's type methods obtaining new derivative-free algorithms. 

Numerical examples that illustrate the theoretical results presented in this paper are also given. We emphasize the 
necessity sometimes of a new way to compute the matrix of divided difference operator since in many examples the 
theoretical one and the classical way used until now to calculate it do not coincide. With the new definition the 
estimate of the computational efficiency index is revisited. Of course, other possible examples of divided differences 
should be considered. 

An approximation of the root with high precision is obtained if we use a multi-precision arithmetics. To illustrate 
the technique presented three examples are worked an completely solved. An approximation of the computational 
order of convergence is computed independently of the knowledge of the root and it confirms in a numerical way 
that the convergence order of the methods is well deduced. 
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